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. ABSTRACT 

Context. Isolated cooling neutron stars with thermal X-ray emission remain rarely detected objects despite many 
searches investigating the ROSAT data. 
. Aims. We simulate the population of close-by young cooling neutron stars to explain the current observational results. 

' Given the inhomogeneity of the neutron star distribution on the sky it is particularly interesting to identify promising 

sky regions with respect to on-going and future searches. 

Methods. Applying a population synthesis model the inhomogeneity of the progenitor distribution and the inhomogeneity 
. of the X-ray absorbing interstellar medium are considered for the first time. The total number of observable neutron 

Oh' stars is derived with respect to ROSAT count rates. In addition, we present sky maps of neutron star locations and 

discuss age and distance distributions of the simulated neutron stars. Implications for future searches are discussed. 
Results. With our advanced model we can successfully explain the observed log N - log S distribution of close-by 
neutron stars. Cooling neutron stars will be most abundant in the directions of rich OB associations. New candidates 
^ ' are expected to be identified behind the Gould Belt, in particular in the Cygnus-Cepheus region. They are expected to 

be on average younger and then hotter than the known population of isolated cooling neutron stars. In addition, we 
propose to use data on runaway stars to search for more radio-quiet cooling neutron stars. 

> 

lO ■ 1- Correction remark - added March 2010 2. Introduction 

^ ■ 

^ • Working on the new version of our population synthesis Mo^g than 10 years after the discovery of its brightest mem- 
^ . code, we recently discovered a bug in the version used in ^er RX J1856-3754 (Walter et al. 1996), the small group of 
■ Posselt et al. (2008). This bug concerns the application of ggven radio-quiet isolated neutron stars (NSs) detected by 
(3 ; the improved ISM models to calculate the absorbed X-ray the ROSAT satellite have gained an important place in the 
, flux. The code erroneously did not cover the whole Galactic ^.j^h zoo of compact objects. Together with Geminga and 
. ^ ' coordinate {l,b) range of the ISM-models to obtain the ab- several close-by young radio pulsars, these objects form the 
sorbing column density N{H). Instead, only a small region ^Q^al population of cooling NSs. Studies of this group of 
^ ( AZ w 7deg, A6 « 7deg ) around the Galactic Center sources already provided a wealth of informat i on on NSs 
; was used, which led to a significant overestimation of the physics (see e.g. baberll [20071: iPaee et al.l [20071: [Z^ [20071 
absorption translating into an underestimation of the pre- \qj. recent reviews) 
dieted neutron star number. 

-in . 1 ■ • Tj J J Since 2001 the number of known close-by radio-quiet 

While all our mam conclusions remain valid we add -.,„ , , , . , ., „ , , -i ,-r 

. J J. -1 J- ii J. J ij. • J.1 iNSs has not been growing despite all attempts to identity 

notes on some details regarding the corrected results m the , . , . ^, Tlr^n Am i ^ \ ^■ ^ ^ r, ^^ ^^ ■ • 

A J- m TT J i J c 1 • 4-u A J- m (mamiy m the RUbAl data) new candidates. Partly this IS 

AppendixLAl Updated ngures are shown m the AppendixLAI \ ^ li n ^ ^ ,i , i i- i ■ ^ 

„ T, . ,• J ° u ur ■ ii, i J- u j-u due to the fact that all these searches are blind. I.e. not nec- 
as well. It IS mentioned m boldface m the text when there ., ^ , , ... . m i , i 

1 , , ., , 1 essarily targeted to promising sky regions, lo advance the 

IS an updated result available. ■ ^ ^-^ ^- r i i- • 

identmcation ot new near-by cooling iNSs it is necessary to 

We note that this is the only paper, where the presented perform a realistic modeling of this population. The seven 

results are influenced by the bug m our population synthesis rqSAT radio-quiet NSs form the main part of the known 

code. Previous works don2t use it, and subsequent works j^^^i population of cooling NSs. For these seven sources 

like, e.g., [ Popov et al. | (|2010D have already applied the cor- ^^^.^ Magnificent Seven. There is still no gen- 

rected version of the population synthesis code. ^^^1 ^^^^ ^^^^^^^ the whole population to which these 

sources belong (terms such as XDINS - X-ray dim isolated 

Send offprint requests to: B. Posseh NSs, DINS - dim isolated NSs, RINS - ROSAT (or ra- 
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dio quiet) isolated NSs, are used sometimes by different 
authors). In this paper we use the term ICoNSs (Isolated 
Cooling NSs). In our opinion, this variant better reflects 
the nature of these sources than, for example, XDINS, as 
many of them are not dim in X-rays. Finally, we call coolers 
all NSs, whose residual thermal emission can be detected. 
Besides the radio-quiet NSs, such sources can be also known 
as normal radio pulsars, or they may demonstrate some 
kind of radio emission or even 7-ray emission distinct from 
classical pulsars. In this study we do not care if an object 
shows some activity in addition to thermal surface emission. 
The only crucial point is the detection of X-ray emission re- 
lated to the cooling of initially hot compact objects. 

The population synthesis (PS) approach is used here 
to investigate the population of close-by young cool- 
ing N Ss. There are two main versions of this method 
(|Fritze-v. Alvenslebedi2000l ). One, which sometimes is also 
called empirical population synthesis, is a kind of top-down 
approach, when objects, for which only integral characteris- 
tics are available, are studied. Numerical models are applied 
to reproduce observed properties of such sources by using 
different sets of subpopulations which constitute an object 
class. For example, the stellar content of a galaxy can be 
derived by modeling its integral spectra. 

Here, we focus on the second - bottom-up - approach, 
often called evolutionary synthesis. In this case the evo- 
lution of a population of objects is modelled from their 
birth based on defined initial distributions of parameters 
and some equations which describe evolution of these dis- 
tributions in time. The main reasons to use this second 
technique are the following. At first, if initial parameters 
or/and evolutionary laws for some kind of sources are not 
well known, then this kind of PS can be applied to test hy- 
potheses about these uncertain properties by comparison 
between modeled and observed populations. Determination 
of the initial parameters of radio pulsars can be a good ex- 
ample here (jFaucher-Giguere fc Kaspill2006[ ). Next, even if 
the initial distributions and their evolution are known, but 
at the moment just a small part of the population - a tip of 
an iceberg - is observed, then PS calculations can be used 
to predict properties of the unobserved faint part of the 
population, and to plan a search strategy to identify new 
members. A more detaile d review of the PS techniqu e can 
be found, for example, in iPopov &: Prokhorovl ()2004f ). 

In this paper an upgraded population synthesis model is 
discussed for the population of close- by (< 3 kpc) isolated 
NSs which can be observed via their thermal emission in 
soft X-rays. Previously, our models were applied to confirm 
the link between th e Magnificent Seven and the Gould Belt 
(jPopov et al.] 20031) and to test theories of thermal evolu- 
tion of NSs ()Popov et al.ll2006[ ). The major interest of the 
present study is to understand how to find more objects of 
this type. After the description of the PS model ingredients 
in Sec. |31 we present our results - log N - log S-curves, age 
and distance diagrams as well as sky maps of expected NS 
locations ~ in Sec. 31 and discuss possible applications in 
Sec. 14.41 An outlook concerning the search for ICoNSs is 
given in Sec. l4.5I Finallv. we give our conclusions in Sec. [5] 



3. The new population synthesis model 

The main physical ingredients which constitute our popu- 
lation synthesis model are the following: 



(A) the initial NS spatial distribution and the NS birth rate; 

(B) the kick velocity distribution of the NSs; 

(C) the Galactic gravitational potential; 

(D) the distribution of NS masses, in the following: NS mass 
spectrum; 

(E) the NS cooling curve; 

(F) the NS surface emission in X-rays; 

(G) the interstellar absorption of X-rays; 

(H) the properties of the X-ray detector. 

At each time step in the simulation we consider at the 
same time eight different NS masses with corresponding 
cooling curves, following each NS on its trajectory through 
space and time until it is too faint to be observed. The 
overall result of, e.g. 50000 x 8 simulated NS evolutionary 
tracks from birth till the time when the temperatures falls 
below 10^ K, is normalised by the mass distribution as well 
as by birth rates (see below for both). Ingredients B, C, E, 
and F are un chang ed wit h resp ect to our previous studies 
(|Popov et all 120031 l2005l l200a hereafter Paper I, II, III), 
so we just briefly comment on them below. The ingredients 
A, D, G, and H - the initial spatial distribution, the NS 
mass spectrum, the interstellar absorption, and the X-ray 
detector properties - are modified, and they are described 
in subsections below in more details. 

For the kick velocity distribution ( B) of newborn neu- 
tron s tars we use that proposed by lArzoumanian et"aI1 
(l2002h . as in all our previous studies. This is a bimodal 
distribution, consisting of two Maxwellians with most 
probable velocities V^caki ~ 127 km and Vpcak2 ~ 
707 km s~^, the average velocity is w 527 km s~^. Recently, 
this assumption was criticized. A single-component dis- 
tributi on with Vueak ~ 400 km s~^ was advocated 
bv ,Hobbs et al.l ([200^, and by iFaucher-Gieuere fc Kaspil 
(|2006[ ) who accepted a non-Maxwellian velocity distribu- 
tion. However, since all objects under study are very young 
and so still are relatively close to their birth place, the ex- 
act shape of the velocity distribution does not have a strong 
infiuence on our final results. 

The calculation of NS trajectories in the Galactic grav- 
itational potential (C) follows that presented in Papers I 
and II . As before we use the axisymmetric Miyamoto-Nagai 
(1975) potential with disc, buldge and halo contributions. 
This potential (with some modifications) was actively used 
in NS calculations by different authors alre ady in the 1990s 
(see, for example, iBlaes fc Madaul (|1993D and references 
therein). Again, it is noted that since the NS are young 
the use of a very precise multicomponent potential is not 
necessary. Thus, we neglect also the Galactic bar which ex- 
tends to a maximum Galactocentric distance of 3.5 kpc (e.g. 
iPichardo et al.|[2004[ ) because its influence is very small for 
the relatively nearby objecs of interest. 

The thermal evolution of NSs is one of the main ingre- 
dients of the study. There are strong uncertainties in this 
topic. Sets of curves have been calculated by many research 
groups for different models of NS interiors (see a reviews, 
e.g. by Page & Reddy 2003 or lYakovlev fc' Pethick 2004|). 
In the following we consider only hadron stars and use cool- 
ing curves (E) corresponding to those marked as model I 
in Paper III . This cooling model includes the following 
important ingredients. Supcrfluid gaps are taken according 
tolTakatsuka & Tamagaki (2004|). Medium modifications of 
the neutrino processes are taken into account. Formation of 
the pion condensate is possible. The relation between tem- 
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peratures of the outermos t core layer and the surface is 
taken from iBlaschke et al.l (|2004l ). No additional heating is 
used. 

In our calculations we do not account for atmospheric 
reprocessing of the thermal radiation, and assume that the 
emitted spectrum is a pure blackbody (see a detailed de- 
scription of the modeling of the observed neutron star ther- 
mal emission in Zavlin 2007). Unfortunately, the distribu- 
tion of the chemical compositions of atmospheres is not 
known for the whole population of young NSs. If the atmo- 
spheric composition changes in time or not is also unclear. 
As an additional parameter one needs further to consider 
the (unknown) magnetic field strength and its effects on 
the atmosphere. So, the available data are not complete 
enough and statistics is too low to consider atmospheres in 
population synthesis studies. 

However, isotropic blackbody emission is a simplifica- 
tion since there is growing evidence that the thermal X- 
ray emission observed in case of the Magnificent Seven 
actually comes f rom hot spo ts and not from the whole 
NS surface (e.g. iHaberl et al.i .200 6; Schwope et al.l 120051: 
iTrtimper et al.l l2004[ ). This is expected to influence not 
only the cooling times but also the detection probabilities. 
But the overall quantitative effect is currently unclear as 
at first it is necessary to include non-uniform heat trans- 
port in the outer NS layers into cooling curves calculations 
self-consistently, for ex ample, in the way it is discussed in 
lAguilera et al.l ((loOT^) ■ Strong toroidal magnetic fields pro- 
duce hot spots and can keep them hot over a longer time 
than th e usual lifetime ex pected from the previous cooling 
models (|Page et al.ll2007D . 

Since the hot spot of a neutron star will have a higher 
temperature than in the case of isotropic cooling, its 
detection would be less affected by interstellar absorp- 
tion. The fact that these effects are not accounted for by 
our applied cooling curves leads to underestimating the 
number of observab l e neu tron stars. On the other hand, 
lYak ovlev fc Pethickl (j2004j) have pointed out that at high 
magnetic field strength the break of the cooling curves 
(marking the transition from neutrino to photon cooling) 
shifts towards lower ages, viz. from ^10^ years at B'-^ OG 
to a few times 10^ years for magnetic fields typical for the 
Magnificent Seven, B-- 10^^ - 10" G. Neglecting this effect 
will lead to an overestimate of the number of observable 
neutron stars with ages of (0.5-1) x 10^ years where some of 
the Magnificent Seven are found . Magnetic field decay and 
the additional heating due to it (jAguilera et al.l[2007bf) can 
further increase the number of observable NSs. We cannot 
make a firm prediction about the net change to the log N 

- logS which is caused by these effects which are partly 
counterbalancing. However, we believe that the blackbody 
estimate of the count rate S is reliable within a factor of 

- 2. 

The soft X-ray spectra of objects under consideration 
are well fitted by the blackbody shape, and, more impor- 
tant, the cooling curves used here were calculated under the 
assumption of blackbody emission of the whole NS surface, 
and cooling curves for non-uniform emission are currently 
not publicly available. We therefore stick to the cooling 
curves of paper HI and the assumption of spherical sym- 
metric blackbody emission, and apply them as a reference 
model to study the main aspect of this paper - the influence 
of the inhomogeneous distributions of neutron star birth 
places and interstellar absorptions in the Galaxy. 



3.1. Initial spatial distribution of NSs 

Following the results of previous investigations (paper I), 
we take as established that the population of nearby NSs 
is genetically related to the Gould Belt. The contribution 
of the Gould Belt dominated the production of compact 
rem nants in the s olar proximity over the past ^ 30 Myrs 
(see iPoppell Il997l for a detailed description of the Gould 
Belt structure). About two thirds of massive stars in the 
^ 600 pc aroun d the Sun belong to the Gould Belt 
(jTorra et al.ll2000D . 

In paper I and II we considered that NSs are born ei- 
ther in the Belt or in the Galactic disc, both treated as in- 
finitesimally thin, disc-like structures, with a uniform spa- 
tial distribution. The total NS birth rate in the Galactic 
disc (the whole area inside 3 kpc) and that in the Gould 
Belt alone were taken to be 270 Myr~^ and 20 Myr~^, re- 
spectivel y. The first value corres ponds to the value inside 
1 kpc by Tarnmann et al.l (jl994 ). The latter value is from 
iGren ier' (20001) ■ and appears in reasonable agreement with 
the hi stori cal supernov a rate inside 1 kpc, as estimated by 
ITammann et al.l (|1994D . 

Here, the previously simplified progenitor distribution 
is significantly upgraded by adopting a more realistic de- 
scription of the massive star geography in the solar neigh- 
borhood. In particular, in the innermost 500 pc from the 
Sun the initial spatial distribution of NSs is assumed to 
follow the presently observed distribution of massive stars. 
Outside the well-known 500 pc and up to 3 kpc, NSs are 
born either in the Galactic disc, assumed to be exponential 
both in z and R directions, or originate from known, rich 
OB associations. 



3.1.1. Inside 500 pc 

For the stars inside 500 pc we use th e Hipparcos data on 
massive stars (jPerrvman fc ESAlll997f) . Our assumption is 
that the present day distribution of massive stars represents 
their distribution in the last 1-2 million years well. 

All the 570 classified B2-08 stars with known parallaxes 
> 0".002 are considered. According to the birthrates above 
in our model 27 (out of 270) NS producing supernova events 
per Myr are assumed to be related to this population. This 
number should be compared with 20 NSs assumed to be 
born in the entire Gould Belt in 1 Myr in our old model. 
The parameters of the new model roughly lie in between 
the old one for i?Bcit = 300 and 500 pc. 

In the code the initial position of each newborn NS com- 
ing out of this population coincides with the position of one 
of these 570 Hipparcos stars. This assumes that the spa- 
tial distribution of the massive main sequence stars does 
not change significantly over the NS evolutionary timescale 
(w 1 Myr). A progenitor star is chosen randomly, but we 
tried to account for the probability for each particular star 
to explode. Since for each star the spectral type is known, 
one can roughly estimate its mass. The mass interval corre- 
sponding to each spectral type gives us an estimated range 
for its total lifetime. If one has no information about the 
age of the star then the probability is assumed to depend 
just on the star's full lifetime: the probability is higher for 
more massive stars as they are shorter lived (actually, the 
probability to explode is inversely proportional to the life- 
time of the star) . 
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Fig. 1. The new NS progenitor distribution. The grey triangles represent the stars in the OB associations. Crosses stand 
for the considered Hipparcos stars inside 500 pc. Smah circles mark random stars in the Galactic disc (500 < d < 3000 pc). 
Outside 500 pc from the Sun two populations are visible: NSs formed in OB associations, and NSs formed in the field. 
Bottom panels show the projected distribution X-zoomed to 500 pc around the Sun. Two circles on the bottom left panel 
are plotted with radii of 200 and 400 pc from the Sun. 



As noted above, the mass interval corresponding to each (ii) 
spectral class is converted into an approximate limiting life- 
time: Tlm is the longest lifetime for the lowest mass, and 
T/fM is the shortest lifetime for the high mass end of the 
particular mass interval. If the star is known to belong to an (iii) 
OB association, we also consider the age of the OB associa- 
tion. The ages of th e OB associations are s hown in Table [1] 
they are taken f rom de Zeeuw et al.l ()1999t ). In several cases (iv) 
Ide Zeeuw et al.l (| 19991 ) give iust one number (not a range) 
as an estimate of an age. Of course, these numbers have 
some uncertainty. To account for this we added an uncer- (v) 
tainty (±1 Myr) in all those cases. By this we have for all 
OB associations in the Table TLq b and TVqb , a lower and (vi) 
upper age limit for the OB association, respectively. 

The probability for a progenitor star to explode (and so 
to contribute to NS production in our model) is estimated 
as follows. 



Conversely, if Tlm < TLqb then the probability is in- 
versely proportional to the longest lifetime, T^m, of the 
star. This probability is also used for stars not belonging 
to OB associations within a 500 pc distance. 
If TifM < TLqb but Tlm > TUqb then the probability 
is inversely proportional to {TUob — TLqb), as one has 
to consider overlapping time scales. 
If Thm < TLqb and TLqb < Tlm < TUqb then 
the probability of exploding is inversely proportional to 
(Tlm — TLqb)- 

If TLqb < Thm < TUob and Tlm > TUob then the 
probability is inversely proportional to {TUob — Thm)- 
If TLoB < Thm, but Tlm < TUob then the proba- 
bility of exploding is inversely proportional to the min- 
imum time value {Tlm — Thm)min = 1-44 Myr, which 
was obtained for the spectral type BO corresponding to 
the mass interval [15 Mq, 18 Mq]. 



(i) If Thm > TUob, then the probability of a star explod- 
ing is zero, as the shortest estimated lifetime is longer 
than the age of the OB association. Such stars do not 
contribute to NS production in our model, i.e. no com- 
pact objects could have been born at their locations 
because the OB association is too young. 



To summarize, the initial spatial distribution of NSs 
inside 500 pc in our model follows the distribution of 
Hipparcos massive stars weighted by their explosion prob- 
abilities as estimated above. Outside 500 pc this weighting 
is not applied. 
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3.1.2. Outside 500 pc, up to 3 kpc 

Most of the NSs in our model, 200 out of 270 per Myr, 
are born in one of 49 OB associations outside 500 pc and 
up to 3 kpc. The radii of the association extents are taken 
to be equal to 100 pc. This is a safe upper limit on their 
sizes , as usually diameters ar e about 80 pc, and up to 200 
pc (Me rnik fc EfremovlfToOSl) . The probability to be born 
inside this radius is uniform. We assume that the number 
of stars with known photometry roughly reflects the total 
number of stars in the cluster. The probability to be born in 
each cluster is proportional to this quantity. The data on 
the 49 OB associations are taken fro m the catalogue by 
Blaha fc Humphreys! (|1989D (see also 'Mel'nik fc Efremovl 
19951) . Note, that distances to the a ssociations are reduc ed 
by 20% according to conclusions bv lDambis et al.l ()200lD . 

Finally, out of the total 270 NSs, 43 are supposed to 
be born during 1 Myr in the Galactic disc outside 500 pc 
from the Sun, but inside 3 kpc. Their positions in the disc 
are calculated randomly. The disc is exponential both in R 
and Z. For the disc Z-scale we take the value 100 pc. For 
the radi al scale we take the value 3 kpc ()Kent et al.lll99ll 
iFux fc M artinet 19941 . 

In our model the NS formation rate outside 500 pc 
around the Sun corresponds to ~ 0.884 10~^^ pc~^ yr~^. 
This r ate is lower than the one stated by iTammann et al.l 
()1994[) . They, adopting that the Milky Way is type Sbc, 
and combining the corresponding value with the histor- 
ical rate in the Galaxy, reported a lo cal SN rate of 
2 10~^^ pc~^ vr~^. iTammann et ahl ()1994l) gave also a rate 
of 2.9 10^^^ pc^^ yr~^ which is based on the census of O- 
B2 stars in 1 kpc around the Sun. Clearly, both values 
should be dominated by the productive Gould Belt. Outside 
the Belt the rate should be lower, however, it should be 
anisotropic and non-uniform on the scale of 3 kpc. Our 
adopted values, ~ 3.438 10~^^ pc~^ yr~^ inside 500 pc and 
~ 0.884 10~^^ pc~^ yr^^ outside 500 pc are chosen to reflect 
this fact. At the 500 pc transition our model gives compara- 
ble values due to the scaling through OB associations. We 
checked that at this transition our model distance distri- 
bution of unabsorbed NS show consistently increasing num- 
bers. Probably, we slightly underestimate the NS formation 
rate at larger distances in the directions of the Galactic 
center and spiral arms. However, for our calculations closer 
regions are more important (see distance distributions of 



OB association 


Age [MyrsJ 


Upper Scorpius 


4-6 


Upper Centaurus Lupus 


12-14 


Lower Centaurus Crux 


9-11 


Vel OB2 


6-10 


Trumpler 10 


15-30 


Collinder 121 


4-6 


Ori OBI 


9.5-13.3 


Per OB2 


3-7 


a Persei (Per OB3) 


25-50 


Cas-Tau 


25-50 


Lac OBI 


2-25 


Cep 0B2 


5-10 


Cap 0B6 


25-50 



Table 1. OB associations inside 500 pc and their ages (from 
de Zeeuw et al. 1999). For Ori OBI the age proposed for 
its older part Ori OBIA is taken. 



ICoNSs below). Anyway, we warn the reader that some un- 
certainty in the birth rate normalization is hidden inside 
our model. Comparison of the Log N - Log S distribution 
for old and new models of the initial spatial distribution of 
NSs is made in Fig. [3] and discussed in SecU) 

3.2. The interstellar absorption of X-rays 

X-rays are absorbed by the heavy elements of the inter- 
stellar medium. This is especially true for the soft X-rays 
below 1 keV. Their extinction has to be accounted for if 
one aims at realistic estimates of the sources observabil- 
ity. The amount of material between the object and the 
observer is measured as hydrogen column density N{H). 
Then, assuming the same abundances towards all direc- 
tions and accounting for photoelectric absorption processes 
by atoms, ions and dust particles one can determine the 
X-ray extinction. The treatment of X-ray absorption is sig- 
nificantly upgraded in our new population synthesis with 
respect to previous studies by considering all these aspects 
, i.e. a more realistic model for N( H), up-to-date abun- 
dance tables and cross sections (see iPosselt et al.l ()2007() 
for a more thorough discussion and paper HI for details 
about previous treatment of the ISM absorptiorQ ) . 

In the following we us e the abundance tab les and pho- 
toelectric cross sections bv I Wilms et al.l (2000) as combined 
in the tbabs routine which is implemented e.g. in XSPEC0. 
We also introduce two different, new N{H) models in three 
dimensions to account for the highly inhomogeneous distri- 
bution of the interstellar medium (ISM) . Both models take 
into account the most upd ated results in th e loca. 1 neigh- 
borhood (up to 230 pc) bv iLallement et al.l ^qM)- Their 
work well describes the local bubble which is surrounded by 
clumps of denser ISM and connected by chimneys to other 
bubbles. For distances larger t han 230 pc w e use t he ana- 
lytical ISM model described in iPopov et all ()2000D m one 
of our new N{H) models, in the following it is named "ana- 
lytical ISM model" . For the second N(H) model we use the 
extinction study bv iHakkila et al.l (|1997 ^. In the following 
we refer to this model as the "Hakkila ISM model" . 

Concerning the transition region at 230 pc we want 
to make the following short comments. The column den- 
sity can only increase with larger distance, thus the new 
data from the analytical formulae / the Hakkila study are 
only applied if the resulting column d ensities are larger 
than those from ILallement et al.l ()2003l) at 230 pc. At dis- 
tances larger than 230 pc the models become naturally 
much coarser due to fewer measurements. For the analytical 
model the transition is smooth since the analytical formulae 
at, e.g., 240 pc resii l t in v alues similar to the ones at 230 pc 
bv lLallement et all (j2003|) . In case of the Hakkila model the 
transition is more abrupt due to the patchiness of the avail- 
able measurement d ata applied for the extinction study by 
IHakkila et all (|1997D . On average the hydrogen colum den- 
sity at around 230 pc is around 10^° cm~^. This average 
value increases by around 10^° cm~^ at the transition for 



We also corrected in the old ISM model the use of full widths 
of half maximum instead of the required v a riances in the expo- 
nential density factors bv' Zane et al.l l|1995f ): [bickev fc LockmanI 
(| 19901 . However, the effect of smaller ISM density above the 
Galactic plane on our model results was quite low. 
^ http:/ /heasarc. gsfc.nasa.gov/docs/xanadu/xspec/ 
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the Hakkila model, while in case of the analytical model 
the average increase is only around 10^^ cm~^. 

Both models were checked with open cluster e xtincti on 
dat a which were presented by iKharchenko et alj ()2005al lbh 
and iPiskunov et alJ ()2006[ ) and are complete up to a dis- 
tance of 850 pc. The data o f the open clusters w ere then in- 
cluded in both models (see 'Posse lt et alJ (|2007t ) for details 
on the models). The column densities derived from both 
ISM mod els at large distances were also compared with the 
results bv lSchlegel et al.l (|1998D from infrared dust emission 
measurements. In both cases, the analytical ISM model and 
the Hakkila ISM model, data cubes in Galactic coordinates 
(steps of one degree) and distance (steps of 10 pc) are used 
for our population synthesis. The sampling of the data cube 
does not reflect the actual accuracy, which is around 25 pc 
within 230 pc and much coarser at larger distances. 

It is noted that the both models are large-scale models 
and even the Hakkila ISM model does not account for small 
ISM clumps, especially at large distances, above 1 kpc. 

3.3. Mass spectrum 

NS cooling curves strongly depend on masses. Once the 
cooling scenario is chosen, it is therefore necessary to spec- 
ify the NS's mass spectrum. As in our previous studies, we 
use the mass spectrum of NS which is derived by the joint 
use of Hipparcos data on spectral classes of clo se-by massive 
stars and calculations bv lWooslev et al.l (|2002[ ). In addition, 
in this paper a modified mass spectrum is s t udied obtained 
by using the calculations from lHeger et al.l (|2005l) . 

Both mass spectra are obtained in a similar way, and a 
priori it is impossible to say that the "new" one is prefer- 
able with respect to the "old" one. We take spectral types 
of massive stars around the Sun as given in the Hipparcos 
catalogue (jPerrvman fc ESA||1997I) . To each spectral type 
we assoc iate a range of masses. Applying now the calcula- 
tions by IWooslev et al.l ([2002) or, instead, bv iHeger et al.l 
(j2005|) we obtain the baryonic masses of compact objects. 
Then the baryonic masses are translated into gravitational 
ma sses by using the simp le equation, valid for hadronic NSs, 
by (jTimmes et al.lll996l) 

Mbar - Afgrav = 0.075 M^.^^. 

The accuracy of this relation is sufficient for our purposes, 
as one needs to know the mass within a few percent (less 
than half the size of a mass bin). Finally, we have a mass 
spectrum with eight mass bins which correspond to the 
eight cooling curves we use in our calculations. 

In Fig. |2| the used two mass spectra are shown - the 
old (dotted line in the graph) mass spectrum which we 
also used in previous studies, and the new one. Note the 
non-equal widths of the mass bins. Eight cooling curves for 
gravitational masses 1.1, 1.25, 1.32, 1.4, 1.48, 1.6, 1.7, and 
1.76 Mq are taken, both in the new and old mass spectrum. 
Each bin boundary corresponds to the mean between two 
values of mass used in calculations. 

For the old spectrum in Fig. |5| we use the relation be- 
tween a proge nitor mass and NS b aryonic mass derived 
from Fig. 17 bv lWooslev et al.l ()2002D . We approximate the 
data in this figure with two linear relations (for progenitor 
masses < 15 Mq and > 25 M©, respectively), and constant 
between 15 and 25 solar masses (see Paper II). A peak at 
1.4 Mq appears due to the fact that in this region pro- 
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Fig. 2. The old (as used in in previous studies, dotted line) 
and the new (as introduced in this work) mass spectra, 
binned over eight intervals of different widths (see the text). 



genitors of different mass produce NSs of similar masses. 
The new mass spectrum, shown with a solid li ne, is derived 
much in the same way, but using data from IHeger et al.l 
(l2005 f) for stars with masses > 12 M©. As it is visible in 
the figure, the peak at 1.4 M© disappears. This is because 
the dependence of the compact star mass on the progenitor 
mass contains no flat regions. The most important feature 
of both mass spectra is a very small number of NSs with 
masses above 1.5 Mq. As the cooling of a NS is faster for 
massive stars, a small number of massive objects implies a 
small number of cold stars in the same age group. 

3.4. Detector properties 

As before (e.g. in Paper I) we use the instrumental response 
of ROSAT to draw conclusions with respect to the ROSAT 
All Sky Survey (RASS). However, the energy range consid- 
ered, e.g. in Paper I was 0.08 keV to 3 keV while the RASS 
data analysis was done for a slightly smaller energy range. 
To be in full agreement we changed the energy range now 
to 0.1 keV to 2 keV. The effective area of the instrument 
is discontinuous over energy due to sharp instrumental ab- 
sorption edges in the energy band of interest. With respect 
to previous studies our sampling of the effective area mea- 
surements was improved to better account for such absorp- 
tion edges. In few cases we found differences of up to 40% 
in the resulting absorbed fluxes for an object. However, this 
effect is pronounced only for very soft sources at very small 
distances. Thus the influence on the overall results of the 
population synthesis are small. We note, that of course one 
has to consider vignetting of the X-rays as not all RASS 
sources are observed on-axis. This is done by an effective 
area measurement with averaged vignetting applicable for 
the RASS. 
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Fig. 3. Log N - Log S for the new mitial spatial distribu- 
tion (solid line) and two variants of the old spatial distri- 
bution with i?bcit — 300 pc (dotted) and 500 pc (dashed). 
The curve shown as a solid line is used as reference on 
the following plots, too. This curve is calculated for the old 
mass spectrum, old ISM distribution, and new ISM element 
abundances. The data points represent the Log N - Log S 
distribution for the known young close-by isolated neutron 
stars with thermal X-ray emission (the Magnificent Seven, 
Geminga, the "second Geminga" - 3EG J1835-I-5918, and 
three radio pulsars: B0833-45, B0656-hl4, B1055-52). Filled 
symbols correspond to the Magnificent Seven, open sym- 
bols to the other NSs. Erro r bars are plo t ted fo r the poisso- 
nian errors. From the list in lPoDOv et all (|2003f ) we removed 
PSR1929+10 since results o f the re cent XMM-Newton ob- 
servations by iBecker et al.l ()2006[ ) demonstrate that the 
spectrum is non-thermal. 



4. Results of calculations and discussion on the 
influence of model modifications 

In this section we compare results of different modifications 
of the population synthesis scenario studied in this paper. 
We start by comparing the log N - log S distributions cal- 
culated with different assumptions, proceed to the resulting 
map of the cooling neutron star distribution on the sky, fol- 
lowed by age and distance distribution, and finally discuss 
the implications on searches for new ICoNS candidates. 

4.1. Comparison of log N - log S distributions for different 
model modifications 

The effects of the successive modifications in our PS-model 
on the log N - log S distribution are shown in Figsl3][S] The 
solid curve in all the figures is the same. It corresponds to 
the new initial spatial distribution of the NSs, new abun- 
dances and cross-sections, the old mass spectrum, and the 
old ISM model. 

At first we discuss the influence of the new initial spatial 
distribution of NSs (i.e. distribution of progenitors). We 
compare three models: two variants of the old model with 
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Fig. 4. In this figure we present the effect of the new mass 
spectrum and new ISM element abundances. All curves are 
plotted for the new initial spatial distribution of NSs. Solid 
curve: old mass spectrum and new abundances. Dotted line: 
old mass spectrum and old abundances. Dashed curve: new 
mass spectrum and new abundances. All curves are calcu- 
lated for the old ISM distribution model. 



RbcU — 300 pc and 500 pc and the model with the new 
initial spatial distribution described in Sec. 13. II Log N - Log 
S curves are shown in Fig. |31 Except for the initial spatial 
distribution all other parameters are the same: old mass 
spectrum (Fig.[2j dotted curve), new abundances, old ISM 
distribution (calculations with the new ISM distribution 
take long times, so for comparison of different modifications 
we use the old ISM distribution model when possible). 

The solid curve corresponds to the new initial spatial 
distribution. The result with the new progenitor distribu- 
tion lies in between those of the two old ones for large and 
moderate fluxes. In fact, our new calculation is only slightly 
above the simpler model with a Gould Belt radius of 500 pc. 
For count rates logs' < — 1 the new initial spatial distribu- 
tion results in roughly 0.1 dex more observable sources. 
The new model of the spatial distribution and the old for 
^Boit = 500 pc are in reasonable correspondence with ob- 
servations at the bright end, where it is supposed that our 
knowledge about observable objects is mostly complete (as 
we see below new models of the ISM distribution shift the 
curve down improving the correspondence). However, the 
locations of observable sources on the sky for the old and 
new initial spatial distributions are significantly different 
(compare Fig. [5] below and Fig. 6 in Paper II ). 

Next, one has to compare the Log N - Log S distri- 
butions calculated for different mass spectra (see Sec. 13. 3p . 
Computations were made for the new initial spatial NS dis- 
tribution, new ISM element abundances and old ISM dis- 
tribution. Results are shown in Fig. U) To see the effect of 
the new mass spectrum one has to compare the solid and 
the dashed curves. The solid curve is the same as in Fig. [3l 
It corresponds to the old mass spectrum (Fig. ^ dotted 
line). The dashed curve, which is nearly indistinguishable 
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from the sohd one, is calculated for the new mass spectrum 
(Fig. [21 solid line). Surprisingly, the two different mass spec- 
tra provide nearly the same result for the Log N - Log S. 
This is very positive for our approach since the mass spec- 
trum is a poorly determined component of the model. As 
noted above, none of the mass spectra applied here is cur- 
rently more preferable than the other. In the following we 
continue to use the "old" mass spectrum. 

In Fig. m the effect of the new ISM element abun- 
dances is also shown. The dotted curve should be com- 
pared with the solid one. As it can be seen, the effect is 
not very important. Naturally, it more strongly influences 
faint sources, which are situated at larger distances (on av- 
erage). However, even for far away sources the effect is not 
significant. 

Finally, we compare the log N - log S curves obtained 
for three different models of the ISM distribution. In Fig. [5] 
we show results for the old simple model (reference curve), 
the new analytical ISM model, and the Hakkila ISM model. 
The curve of the new analytical ISM model is very close to 
the one of the old ISM model for bright sources. As sources 
become fainter the number of observable objects becomes 
smaller than before due to the appearance of directions with 
high absorption. Using the Hakkila ISM model we obtain a 
log N - log S curve which fits the known observational con- 
straints even better. For small fluxes this curve approaches 
the reference one. In comparison to the analytical ISM cube 
there are large "windows" of low extinction in the Hakkila 
ISM model resulting in more sources at low fluxes. The an- 
alytical model on the other hand seems to overpre dict the 
absorption when going to larger distances (Posselt et al.l 
|2007() . Thus, at fainter fluxes one may regard the results of 
the analytical and Hakkila ISM model as lower and upper 
limits on the number of observable sources. 

We conclude, that in terms of the log N - log S distri- 
bution our new, more realistic, models are in good agree- 
ment with the older ones. The interstellar absorption by the 
Hakkila ISM model results in a log N - log S curve situated 
slightly closer to the known observations of INSs. 

By using a different set of cooling curves one may 
improve the agreement between computed and observed 
log N - log S distributions. A detailed analysis is, however, 
outside the scope of the present investigation. Therefore, in 
the following only few comments are made regarding the set 
of cooling curves applied in our previous Paper HI. As ap- 
parent from Fig. [3] and Fig. [5] our final log N - log S curves 
which include the new ISM distribution models are not sig- 
nificantly different from the curve using -Rbcit = 500 pc in 
our old calculations (Paper HI). So, the main conclusions 
about different sets of cooling curves presented in Paper HI 
should not change. While the curve calculated in the case 
of the Hakkila ISM model ~ dashed in Fig. [S]- lies at low 
count rates slightly below the old curve with i?bcit — 500, 
this deviation is not strong enough to validate the cooling 
models which have been rejected in Paper HI. Moreover, 
for this choice of the ISM distribution the cooling models 
which have been considered as acceptable in Paper HI fit 
even better the observations. 

4.2. Sky maps 

In earlier papers on population synthesis calculations we 
did not attempt to produce realistic sky maps of the ICoNS 
distribution, as some of our assumptions were quite crude. 
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Fig. 5. Log N - Log S in case of different X-ray absorb- 
ing ISM models. All curves are plotted for the new initial 
spatial distribution, the old mass spectrum and new ISM el- 
ement abundances. Solid curve: old, simple analytical ISM 
model as e.g. in Paper HI . Dotted curve: new improved 
analytical ISM model; dashed curve: Hakkila ISM model, 
see Sec. 13.21 For an updated version of this figure see 
Fig. lA.ll in the Appendix lAl 



Now, with detailed distributions of ISM and birth places of 
NSs, one can try to compute such maps. In Fig.|B]we present 
a map in Galactic coordinates for the new initial spatial dis- 
tribution of NSs, old mass spectrum and the new analyti- 
cal ISM model. The plot shows the expected number den- 
sity for sources with count rates >0.05 cts s~^. Obviously, 
sources appear to be restricted towards the Galactic plane 
and the plane of the Gould Belt. Few objects are expected 
to be found at latitudes higher than 30° (see, however. 
Sec. 14. 5p . Inside ± 30 degrees from the Galactic plane the 
distribution of sources is dominated by NSs from relatively 
close, rich OB associations. We mark by boxes the most 
important of them: Sco 0B2, Cyg 0B7, Cep 0B3, and Ori 
OBI. Interplay between source distribution and 3D ISM 
structure allows us to make predictions on which directions 
are most promising for looking for new ICoNSs. 

The computed distribution of sources in Fig. [S] is 
dominated by sources with count rates in the range 
[0.05, 0.1] cts s~^. This corresponds to the dimmest sources 
among identified radio and 7-ray silent ICoNS, or even 
to count rates smaller than < 0.1 cts s'-^ where there 
aren't currently any ICoNS known . In the distribution of 
known ICoNSs (jMotch et al.l I2007D there are no sources 
with Galactic longitude in the range ~ 50° to ^ 200°. At 
the first glance, this is in contradiction with the our present 
map. However, it is necessary to note, that for brighter 
sources (> 0.5 cts s~^) the Cygnus-Cepheus region does 
not give a strong contribution , and 1-2 objects could es- 
cape identification in this overcrowded area close to the 
Galactic plane (see Sec. l4.4l for a more detailed discussion). 

In Fig. [7] we present a similar map for the same pa- 
rameters concentrated on the faint sources having nominal 
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Fig. 6. The expected number density of isolated neutron stars with thermal X-ray emission in units of numbers per 
square degree. The Galactic map is in Mollweide projection. Only sources with ROSAT PSPC count rates larger than 
0.05 cts are considered, the same value as used in Paper II in their Fig. 6. The simulation was done for new initial 
progenitor distribution, new abundances, old mass spectrum and the new analytical ISM model, thus corresponds to 
the dotted log N - log S curve in Fig. [S) Marked in blue are the positions of the Magnificent Seven and in black the 
positions of close young radio pulsars with detected thermal X-ray emission. For an updated version of this figure 
see Fig. IA.2I in the Appendix [Al 



ROSAT count rates between 0.001 cts s~^ and 0.01 cts s~^ 
We note that the RASS is not so deep. Fig. [7] demonstrates 
the effect of more distant (> 1 kpc) OB associations, the 
most important ones are marked in the map. 

4.3. Age and distance distributions 

In this paper we introduce age and distance distributions 
of close-by NSs as predicted by our PS model. These distri- 
butions can help to illustrate better the properties of both 
discovered and still elusive ICoNSs. Figures E] and O show 
the computed age and distance distributions for NSs with 
ages up to 10^ years with bin size 10^ years, and distances 
up to 700 pc with bin size of 140 pc, calculated for dif- 
ferent model assumptions ( the two variants of the initial 
spatial distribution of NSs, and the three models of ISM 
distributions) . 

We start by presenting only the general picture and refer 
to Sect. l4.3TT] for a comparison of the impacts the different 
PS scenarios have on our results. Fig. |S] confirms the natu- 
ral expectation that the youngest neutron stars contribute 
most to the observable object fraction in all flux ranges as 
they are both brighter and hotter. The radial distribution 
for different count rates is presented in Fig. |9l While for 
large fluxes most sources are situated in the region ~ 200 - 
400 pc (ie. in the Gould Belt), for fainter objects the picture 
is different. At fluxes below ~ 0.1 cts s~^ one expects to 
see mostly sources behind the Gould Belt. These general 
features are very important, as they indicate that new, still 



unidentified ICoNS are expected to be young objects be- 
hind the Gould Belt. We discuss this proposition below in 
Sec. 1131 

4.3.1. Comparison of the distributions for the different model 
modifications 

The population synthesis models used to obtain the results 
presented in Fig. [S] and Fig. [S] differ on the one hand in 
the initial progenitor distribution and on the other hand 
in the adapted ISM model. As it is already evident from 
the log N - log S curve in Fig. [3] the predicted number of 
observable NSs is generally higher in the case of the new 
initial progenitor distribution compared (white bars) to the 
old one (vertically striped bars). The decrease in number 
from the youngest to the oldest sources for both progen- 
itor distributions is similar. For sources with count rate 
< 1 cts s ~^ the new initial spatial distribution results in 
larger number of young sources. These are objects behind 
the Gould Belt in OB associations. Outside the Gould Belt 
the new distribution is more compact than the old one: 
most of associations are closer than 2 kpc (in the old model 
NSs were uniformly distributed in the disc up to 3 kpc). 

As the next step we compare the results for the old sim- 
ple ISM model (white bars in Fig. |5] and Fig. [3]) with those 
obtained using the new analytical model (black bars). In 
both cases we use the new inital progenitor distribution. In 
general, calculations with the old ISM model produce more 
observable sources (which can also be seen in the log N- 
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Fig. 7. Same as Fig. [51 but here only faint sources with ROSAT PSPC count rates between 0.001 cts s~^ and 0.01 cts 
s~^ are considered. Marked are regions of the corresponding OB associations giving birth to the neutron stars. These 
OB associations have distances between 1 kpc and 2 kpc. For an updated version of this figure see Fig. IA.3I in 
the Appendix lAl 



log S curves in Fig.[S]). The following detailed compar- 
ison of the analytical and Hakkila ISM model for 
specific directions is obsolete after correcting for 
the bug in our population synthesis code and can 
be neglected. The general picture however remains 
unchanged. This difference is more clearly pronounced for 
faint sources with count rates below 1 cts s ~^ and ages be- 
low 3 X 10^ years. The X-ray emission of these objects must 
be more absorbed in the new analytical ISM model to ex- 
plain this effect. Therefore, these young faint sources have 
to be located on directions which correspond to larger hy- 
drogen column densities in the new analytical ISM model. 
Indeed, there are several small regions of higher column 
density in the Galactic plane, and some of them hinder the 
view to star associations, e.g. at / = 215°. 

The distance histograms at different count rates at low 
distances (e.g. up to the third distance bin, [280,420] pc) 
show no difference or only slight reduction comparing the 
resulting numbers from applying the old ISM with those ob- 
tained by using the new analytical ISM model. At larger dis- 
tances in the distance histograms there are fewer observable 
NSs in the case of the new analytical model in comparison 
with the old one, e.g. in the last distance bin, [560, 700 [ pc, 
in case of the faintest sources. These differences are evident 
from the properties of the new analytical ISM model. First, 
due to the inhomogeneities in the new analytical ISM model 
there are regions of higher column density compared to the 
old model. Second, the "average" column density (without 
apparent ISM clumps) in the Galactic plane is less for the 
new model at distances below 450 pc. 

The lower average N(H) value at relatively close dis- 
tances can explain the comparable or even slightly increased 



predicted number of observable neutron stars even with 
present ISM clumps. The ISM clumps are however respon- 
sable for the number reduction in the right part of the dis- 
tance histograms. For example there is strong absorption 
towards the direction I = 100°, b = —5° relatively close to 
the Cygnus-Cepheus OB-associations region. The predicted 
number of observable NSs is reduced there. Since the OB 
associations towards this direction of sight are at distances 
larger than 500 pc (e.g. Cyg 0B7 at 740 pc, Cep 0B2 at 
700 pc and Cep 0B3 at 900 pc. Cam OBI at 900 pc) the 
effect is pronounced mostly at faint count rates (objects 
relatively far away). 

Now we compare the predicted numbers calculated us- 
ing the analytical ISM model (black bars) and the Hakkila 
extinction model (diagonal stripes). Here, we note again 
that up to a distance of 230 pc both ISM models are 
almost the same. F rom 240 pc the extinction study by 
iHakkila et al.l (|l997t) is applied, and so starting from this 
point the absorption in this more inhomogeneous model 
is usually much larger than in the case of the analytical 
model which is characterized by slower and more homo- 
geneous increase of the column density. Only at distances 
around 450 pc do the hydrogen column densities of the an- 
alytical ISM model become larger than those of the Hakkila 
ISM model in most regions. For example, the Hakkila ISM 
model has in general lower hydrogen column densities at 
\b\ > 10° for distances larger than 400 pc. Extinction is 
also low in this ISM model from I = 230° to / = 285° in 
the Galactic plane (all distances), where differences of up 
to an order of magnitude can be reached in comparison to 
the analytical ISM model. 

In agreement with the log N-log S curves of Fig. [S] the 
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Fig. 8. In this figure the age histograms of observabie INSs with thermai X-ray emission are plotted. Three different 
ROSAT count rate intervals represent very bright, bright and moderatly brigh t ROSAT X-ray sources (note the RASS 
Bright Source Catalogue limiting count rate of 0.05 cts s~^: IVoees et al.l [199911 . Different bars correspond to different 
variants of our model. The bars with vertical stripes indicate the results of the old population synthesis model with a 
Gould Belt radius of 500 pc (Paper II , see also Sec. 13. ip . the white bar correponds to the population synthesis applying 
the new inital spatial distribution (Sec. 13.11 reference log N - log S curve in Fig. [3] to Fig. [SJ. Results obtained in the 
frame of the new analytical ISM model are shown with black bars, and those obtained with the Hakkila ISM model are 
represented by bars with diagonal stripes (see Sec. l3.2l and Fig. [5] for the corresponding log N ~ log S curves). Please note 
different y-scales in different panels. For an updated version of this figure see Fig. IA.4I in the Appendix O 



Hakkila NS number of bright sources is slightly less than 
in the case of the analytical model, both in the age and 
the distance histograms. This slightly lower observable NS 
number is due to larger absorption at closer distances. For 
count rates below 1 cts s ~^ the application of the Hakkila 
ISM leads to a higher number of observable NSs, especially 
of young age. The lower absorption at large distances ex- 
plains these higher observable number of ICoNSs. For the 
discussed distance range, however, the distance histograms 
indicate only small differences for sources fainter than 1 cts 
s -\ 



4.4. Resultant search strategies 

The main aim of this study is to make advances in the 
strategy for searching for new ICoNSs. Based on the re- 
sults presented above we can suggest some modifications 
in the approach to look for new candidates. According 
to what we discussed in the previous sections, new can- 
didates expected to be identified at ROSAT count rates 
< 0.1 cts s~^ should be young objects born in rich 
OB associa tions behind the Go u ld Belt. Most of the re - 
ccnt studies 'Agiieros et all (|2006[ ): IChieregato et"all (|2005[ ): 
fRutle dgc et al. (2003j) looked for new candidates far from 
the Galactic plane. It seems that this is not very promis- 
ing. Our results indicate that ICoNS should be searched in 
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Fig. 9. The distance histograms of the observable INSs 
with thermal X-ray emission are shown for three different 
ROSAT count rate intervals. Different bars correspond to 
different variants of our model, and have the same pattern 
as in Fig. [S] Please note different y-scales in different pan- 
els. For an updated version of this figure see Fig. IA.5l 
in the Appendix [Al 



directions of OB associations such as Cyg 0B7 and Cep 
0B3. On average, new candidates should be slightly hotter 
than the known ICoNSs, as they are younger. 

T he absence of sour ces in the second Galactic quad- 
rant (|Motch et al.l [2007h can be naturally explained tak- 
ing into account the properties of close-by OB associa- 
tions. Most of the associations (Per 0B2, Cas-Tau, a Per, 
Cep 0B6) marked in this empty region in Fig. 3 from 
iMotch et al.l (|2007[ ) have ages which do not favor the ap- 
pearance of young observable NSs at present times. They 
are too young (< 7 Myrs) or too old (> 25 Myrs), see 
Table 1 above. As noted above, the strong contribution of 
the Cygnus-Cepheus region to the expected neutron star 
number, shown in map Fig.[BJ comes especially at low count 
rates. Thus, we expect that sources with 50° < 1 < 200° 
would be identified having lower fluxes than the known 
sources, since they come from OB associations further away. 

Considering sky coverage the ROSAT All Sky Survey 
(RASS) is currently the best choice to look for new ICoNSs 
in the Cygnus-Cepheus region which is, according to our 
presented results, the most promising region. However, the 



relatively large positional error circle of ROSAT usually 
includes many possible optical counterparts, especially at 
these low Galactic latitudes. Furthermore one has to ex- 
clude variable X-ray sources to find ICoNSs. In this re- 
spect the recently published XMM-iVewton Slew Survey 
(EPIC with medium filter, e.g. lEsauei et al.ir2006[) may be- 
come an important database. Currently this survey (0.2- 
12 keV band) covers roughly 15 % of the sky, and be- 
low 2 keV its sensitivity is comparable to that of the 
RASS (jFrevberg et al.l [20061 ) . The sensitivity is however 
strongly inhomogeneous perpendicular to the slew direc- 
tion. The astrometric accuracy is around 17" (90% con- 
fidence level), but in few cases an error of 1' due to an 
attitude problem has been reportetj]. In conjunction with 
the RASS, the XMM-Newton Slew Survey with its in- 
creasing coverage can be used to identify non-variable X- 
ray sources and eventually improve their positional accu- 
racy. With the planned expanded energy range down to 
150 eV (M.Freyberg, pers. comm.) the soft X-rays emitted 
by ICoNSs could be detected. However, the yJslWl-N ewton 
Slew Survey has strongly inhomogeneous sensitivity, and 
the RASS is not deep enough for the expected new faint 
sources. Positions of both can be too inaccurate in the 
highly populated Galactic plane. Pointed observations by 
IQAyi-N ewton with its high sensitivity and by Chandra 
with its superb positional accuracy may solve these prob- 
lems for particular directions. 

4.5. Outlook 

Whether the calculated log N-log S curves are realis- 
tic descriptions of the ICoNS numbers also at low count 
rates could be best tested with the help of an All Sky 
Survey like that with the planned cROSITAH. This ex- 
periment is expected to have an angular resolution of 
< 25" and to be 30 times more sensitive than the RASS 
(iPredehl et al.l[2006l) . eROSITA's CCD technology is based 
on an improved concept of the successful XMM-iVewton 
EPIC-pn CCDs. While the noise at low X-ray energies 
is expe cted to be much lowe r compared to the EPIC-pn 
CCDs llMeidinger et al.l [20061 ) the influence of the planned 
eROSITA filter actually reduces the sensitivity to or be- 
low that of ROSAT in case of energies lower than 0.2 keV 
(P. Predehl, pers. comm.). Using a preliminary eROSITA 
response matrix one can estimate that detections of soft 
blackbody spectra are possible with approximately one 
eighth of the RASS count rate after the planned 4 years 
of survey observations. Thus, in principle, sources could be 
expected to reach count rates down to log S« —2.7 (using 
0.017 counts s~^ as the RASS limit). This can result in a sig- 
nificant increase in the number of currently known ICoNSs, 
up to an order of magnitude. However, the actual number of 
identified NSs will depend strongly on the finally reached 
positional accuracy, especially in the Galactic plane. For 
finding new sources of the Magnificent Seven type it will 
be advantageous to have X-ray positions with 1" ac- 
curacy to avoid confusion with faint extragalactic objects 
(|TriimDej|2006[) . 

Besides the standard approach based on X-ray/optical 
wavelengths we mention below other possibilities to look 
for ICoNSs. They are not directly related to our PS results. 



^ http: / / xmm.esac.esa.int / external/xmm_science / slew_survey /xmmsl 
^ www.mpe.mpg.de/erosita/MDD-6.pdf 
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Keeping the restrictions of our model in mind, however, 
they may be important in terms of finding new ICoNSs. 

Recently, ICrawford et al.l (pOOfil found 56 well- 
determined EGRET error boxes devoided of radio pulsars. 
Estimates of the number of young NSs which can man- 
ifest themselves as 7-ray sources show that many of the 
56 sources can be related to these compact objects, but a 
significant frac tion of them should not appear as standar d 
radio pulsars (jHarding et al.l 120071: iGonthier et al.ll2007D . 
Depending on the 7-ray emission model the fraction of 
Geminga-like objects can reach between 3 % a nd 90 % of 
7-ray pulsars according to iHarding et al.l (|2007[ ). They can 
be detectable in soft X-rays as coolers, elusive in the radio 
bancO- New GLAST and AGILE observations can provide 
much smaller error boxes, so soon it could be possible to 
identify sources if they are ICoNSs. It would be fascinat- 
ing to find more "7-ray selected" coolers. Of course, it is 
important to make a joint population synthesis for thermal 
and non-thermal X-ray sources, but this task is out of the 
scope of this paper. 

Another possibility to find new ICoNSs is to search 
for (un)bound compact companions of OB runaway stars. 
More than one hundred OB ru naway stars are know n in 
a 1 kpc region around the Sun (|de Zeeuw et al.l[l999D . In 
comparison with typical OB stars they are characterized by 
large spatial velocities or/and by large distances from the 
Galactic plane. Two main origins of these large velocities 
are currently discussed: dynamical interaction and ex pul- 
sion of a companion in a close binary system (e.g. Blaauwl 
119611 ). The latter case is interesting for the discussion of 
search for new close-by cooling NSs. 

A binary can survive after the first SN explosion in , 
roughly, 10-20 % of cases fe.g. iPopov fc Prokhorovll2006D . 
Then a runaway system consisting of an OB star and a 
compact object (most probably a NS) is expect e d. A y oung 
NS can appear as a radio pulsar. I Saver et all (|l996l ) and 
iPhilp et aTl (|1996[ ) searched for radio pulsar companions of 
~ 40 runaway OB stars. Nothing was found. This result is 
consistent with the assumption that in less than 20% cases 
OB stars have radio pulsar companions. Still, it is inter- 
esting to speculate that runaway massive stars can have 
ICoNS as companions. Then a companion can be identified 
as a source of additional X-ray emission if a NS is younger 
than ^ 10® yrs. It is expected to have a thermal compo- 
nent with T - 50 - 200 eV and L - lO^o - 10^^ erg s~^ 
in m ost standard scenarios of cooling (sec 'Blaschk e et al.l 
(|2004) : Page & Rcddv (2006); Page et al. (2007) and ref- 
erences therein). Our simple estimates show that a cooling 
NS in a binary with an OB star can appear in the propeller 
stage. In this case additional energy output non-related to 
the surface thermal emission can be expected. 

As mentioned above, in our calculations very few NSs 
can reach high Galactic latitudes. However, in our model we 
do not take into account that progenitors themselves can 
be high velocit y objects (ruiiaway stars or so-called hyper- 
velocity stars. iBrown et al.ll2005f) . In this case a NS can 
be born far away from the Galactic plane. This possibil- 



^ GeminRa itsel f was claimed to be detected in radio 
iMalofeev &i Malovl (|1997fl . but the "second Geminga" - 3EG 
J1835-I-5918 - is observed only in 7 and in soft X-rays. 



ity is important i n connection with the recently discovered 
source, Calvera (iRutled ge et al.l[2007D . which could have 
had a high velocity progenitor. Recent strict boun ds on 
the radio emission from Calvera (jHessels et al.l[2007l ) make 
it probable that this source is similar to the Magnificent 
Seven. 

Hyper- velocity stars (iBrown et al.l l2007l) are expected 
to originate mainly from the Galactic center . They acquire 
large velocities after binary disruption in the field of the 
supermassive Galactic black hole. Another possibility is th e 
ejection from young star clusters (jCvaramadze et al]|2007[ ). 
Initial velocities are about 1000 - 3000 km s~^. In 25 Myr 
(maximum lifetime of a NS progenitor) a hyper-velocity 
star can travel up to 50 kpc, more than enough to explain 
Calvera. However, the number of such stars is not high. The 
total rate of production of hyper- velocity stars of all masses 
is between 10"^ and 10~® yr^^ (Brow n et a"Lll20Q7[ ). If we 
take this rate as ^ 10~® yr^^ for NS progenitors then one 
can expect a few ICoNS on very high Galactic latitudes. 
In this case it is not very probable to find one Calvera in 
^ 10 kpc from us. Still, this possibility cannot be excluded. 
The high velocity tail of runaway stars corresponds to 200- 
300 km s~^. In 25 Myr a NS progenitor can reach z ~ 5 kpc 
above the Galactic plane. Potentially, it is enough to explain 
Calvera, but detailed modelling is desirable. 

Of course population synthesis models in general, and 
the one we use in particular, cannot predict new types of 
sources, since our results do not include populations of cool- 
ers absolutely different from those assumed in the model. 
For example, Calvera can be an evolved version of the Gas 
A central co mpact object (CCO , for more information on 
CCOs see e.g. lKaspi et al.ll2"006D . Our model does currently 
not include objects like this (small emitting area, etc.). In 
general, the differences in properties between CCOs and 
close-by coolers indicate that we currently do not under- 
stand initial properties (and perhaps, evolution) of NSs well 
enough. 

A joint population synthesis of all known types of NSs 
would be useful. The presented population synthesis model 
cannot take into account all interesting possibilities related 
to evolution of cooling NSs. So, as usual, surprises are not 
excluded. 

5. Conclusions 

In this paper we presented our new more advanced model 
for population synthesis of close-by cooling NSs. The two, 
sligthly different, mass distributions we consider result in 
nearly the same observable number of NSs. Detailed treat- 
ment of the initial spatial distribution of NS progenitors 
and a detailed ISM structure up to 3 kpc allow us to discuss 
the strategy to look for new ICoNS. Our main results in this 
respect are the following: new candidates are expected to be 
identified behind the Gould Belt, in directions to rich OB 
associations, in particular in the Cygnus-Cepheus region; 
new candidates, on average, are expected to be hotter than 
the known population of ICoNS. Besides the usual approach 
(looking for soft X-ray sources), the search in 'empty' 7-ray 
error boxes or among run-away OB stars may yield new X- 
ray thermally emitting NS candidates. 
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Appendix A: Erratum 

We report in the following on some details regarding the 
corrected results obtained by the upgraded version of our 
population synthesis code. Ah our main conclusions pre- 
sented above remain vahd. 



B. Posselt et al.: Where to look for new ICoNSs 



15 



A.l. Comparison of log N - log S distributions for different 
model modifications 
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Fig. A.l. Log N - Log S for different X-ray absorbing ISM 
models. This figure is the corrected version of Fig. [S] All 
curves are plotted for the new initial spatial distribution, 
the old mass spectrum and new ISM element abundances. 
Solid curve: old, simple analytical ISM model as, e.g., in 
Paper III . Dotted curve: new improved analytical ISM 
model (corrected); dashed curve: Hakkila ISM model (cor- 
rected) . 



2.1), we adopted a birth rate of 27 Myr ^ up to a dis- 
tance of 500 pc (|Grenierj[2000t iTammann et al.lll994 f) and 
of 270 Myr~^ for the who le distance range from to 3000 
pc (jTammann et al.lll994f) . The supernova rate ran ges from 
17 M yr~^ to 27 Myr"^ in the entire Gould Belt i)Grenieij 
|2004[ ). Only 75% to 87% of the core-co llapse supernova pro- 
duce neutron stars ([Heeer et al.ll2003h . Thus, the birth rate 
in the Gould Belt we chose can be overestimated. However, 
a factor of two in the local birth rate uncertainy seems to 
be unlikely. 

As mentioned in the original paper, the main conclu- 
sions about different sets of cooling curves presented in 
Paper III do not change. The chosen set of cooling curves 
in this paper actually represents the best choice of the set 
from Paper III. All other cooling curves from Paper III 
would result in an even higher log N. It is beyond the scope 
of this erratum to identify a new cooling curve set that 
could lead to a log N - log S curve in better agreement 
with observational points. 

We further note the possible effect of statistical fluctua- 
tions on the uncertainty of the observed log N - log S curve, 
since the overall number of (young enough) neutron stars 
in the Gould Belt is small. We will evaluate this effect in 
more detail in a future paper. 

A. 2. Sky maps 

For completeness we show in Figs. IA.2I and IA.3I the cor- 
rected versions of Figs. 6 and 7 from the original paper. The 
general picture remains the same, but on average more NSs 
are expected in accordance to the log N - log S curve. For 
bright sources, a subgroup of Sco 0B2, Lower Centaurus 
Crux (|de Zeeuw et al.lll999f ). is predicted to have more ob- 
servable NSs compared to other Sco OB groups, which was 
not visible in the map obtained with the original version of 
the PS code. 



The log N - log S curves in Fig. 3 and Fig. 4 of the 
original paper remain unchanged. The log N - log S curves 
for the new analytical as well as for the Hakkila ISM model 
in Fig. 5 of the original paper are updated in Fig. lA.ll 
The corrected curves for both ISM models are now sit- 
uated « 0.3 dex above the observational points. The dif- 
ferences between the Hakkila ISM and the improved ISM 
analytical model are smaller than those obtained by the 
original PS-code. However, as before, the application of the 
Hakkila ISM model results in lower N for high count rates 
than obtained by using the analytical ISM model, and in 
higher N for low count rates. Comparison of our new re- 
sults with observations of bright, cooling NSs indicates that 
the model overpredicts the number of NSs by roughly a 
factor of two for both ISM models. The possible reasons 
for this discrepancy are an inadequate treatment of the 
NS birth rate or of their thermal evolution, or other, not 
yet in this paper investigated properties like atmospheres, 
magnetic fields, or statistical fluctuations. Birth rates of 
neutron stars are highly uncertain, espec ially at larger dis- 
tances (see, e.g., the rece nt discussions bv lKeane &: Kranieii 
(|2008l ) andlOfey (|2009[ )V However, since probably nearly 
all the observed XTINSs originate from the Gould Belt we 
have to discuss the local birth rate in the frame of the dif- 
ferences between our PS model predictions and the obser- 
vational measurements. As discussed in the paper (Sect. 



A. 3. Age and distance distributions 

Age and distance diagrams obtained with the corrected 
PS-code show the same main features as those obtained 
with the old code except that the NS numbers of the 
analytical and Hakkila ISM model (the last black and 
diagonal-striped bars in Figs. 8 and 9 of the original paper) 
are larger than before, as is expected from the log N - log S 
curves. 
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Fig. A. 2. This map is the corrected version of Fig. 6, showing the expected number density of isolated neutron stars 
with thermal X-ray emission in units of numbers per square degree. The Galactic map is in Mollweide projection. Only 
sources with ROSAT PSPC count rates larger than 0.05 cts s~^ are considered. Marked in blue are the positions of the 
Magnificent Seven and in black the positions of close young radio pulsars with detected thermal X-ray emission. Note 
that the runs for Figs. 6 and 7 on the one hand, and Figs. IA.2] and lA.3l on the other hand, are separate realisations of our 
code (e.g., radomised birth places); therefore, individual NS "tracks" are not the same. 
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Fig. A. 4. This figure is the corrected version of Fig. [51 It shows the age histograms of observable INSs with thermal 
X-ray emission. Three different ROSAT count rate intervals represent very bright, bright a nd moderatly bright ROSAT 
X-ray sources (note the RASS Bright Source Catalogue limiting count rate of 0.05 cts s~^: IVoges et al.llT999. ). Different 
bars correspond to different variants of our model. The same pattern as in Fig. |8] are used: The bars with vertical 
stripes indicate the results of the old population synthesis model with a Gould Belt radius of 500 pc (Paper II , see also 
Sec. 13. ip . the white bar correponds to the population synthesis applying the new inital spatial distribution (Sec. 13.11 
reference log N - log S curve in Fig. [3] to Fig. [5]) . Results obtained in the frame of the new analytical ISM model 
(corrected) are shown with black bars, and those obtained with the Hakkila ISM model (corrected) are represented by 
bars with diagonal stripes (see Sec. 13.21 and Fig. lA.ll for the corresponding log N - log S curves). Please note different 
y-scales in different panels. 
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Fig. A. 5. This is the corrected version of Fig. ^ the dis- 
tance histograms of the observable INSs with thermal X- 
ray emission for three different ROSAT count rate intervals. 
Different bars correspond to different variants of our model, 
and have the same pattern as in Fig.[S] Please note different 
y-scales in different panels. 



